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1. Selecting the Sampling Period of the LST System 

The objective of this investigation is to determine the effect of 
varying the sampling period on the dynamic response of the sampled- 
data LST system. 

A range of sampling periods is recommended based on the criterion 
that self-sustained oscillations are to be avoided in the LST system. 

The step responses of the LST system are then investigated when various 
sampling periods are used. 

Detailed description of the LST system with the CMG nonlinearity 
is available in the Final Report, CONTINUOUS AND DISCRETE DESCRIBING 
FUNCTION ANALYSIS OF THE LST SYSTEM, January 1, 1974, prepared by the 
authors for NASA, Huntsville, under contract NAS8-29853. In that report 
describing function analyses are applied to the continuous-data and the 
sampled-data models of the LST system with the CMG nonlinearity. It is 
shown that the 9th-order LST system can be closely approximated by a 
4th-order system. 

Two sets of system parameters (System 1 and System 2) were considered 
in the Final Report. The study included in this report is concerned only 
with System 1. 
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It has been established that for System 1, and with y = 1.38 x 10 
for the CMG nonlinearity, self-sustained oscillations will occur if the 
sampling period T exceeds 0.25 sec approximately. 

In order to carry out the discrete describing function analysis for 
the sampled-data system, a sampler is inserted in the nonlinear loop, and 
thus a two-sampler system results. Computer simulation results show that 
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the two-sampler model gives very good predictions on the occurance of self- 
sustained oscillations in the one-sampler system by the discrete describing 
function analysis. 

Figures 1-1 and 1-2 show the zero-input responses of the LST system 

with two samplers; Figures 1-3 and 1-4 show the responses when there is 

only one sampler in the system. In both cases, the initial value of 0y 

(vehicle position) is 5 x 10" rad while all other initial conditions are 

zero. The sampling period is 0.25 sec. With this large sampling period, 

the system actually settles into a self-sustained oscillation with a 

peak value of 0 y approximately equal to 10"^® rad, although this amplitude 

is not visible from the curves of Figures 1-1 through 1-4. As mentioned 

earlier, the sampling period of T = 0.25 sec can be considered as a 

boundary case between stability and instability. 

It is of interest to investigate the step response of the LST system 

-8 

with and without sampling. A step input of amplitude 5 x 10 K Q is applied 

-8 

which yields a final value of 5 * 10 rad for 0y. Figures 1-5 and 1-6 
illustrate the step responses of the continuous-data LST system. As 
expected, the continuous-data LST system with the designated controller 
parameter K q and has a fairly good step response, since it was demonstrated 
that the system has a relative damping ratio of 70% approximately. 

Figures 1-7 through 1-12 show the step responses of the sampled-data 
system with one sampler when T = 0.05, 0.1, and 0.25 sec, respectively. 

Figures 1-13 and 1-14 show the step responses of the two-sampler system 
with T = 0.25 sec! 

When T = 0.25 sec, the step responses again have small oscillations in 
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the steady state. 

The step responses in Figures 1-5 through 1-14 show that the LST 
system with a step input behaves very similar (except for a shift in 
the reference of 8y) to the system with zero input and nonzero initial 
value for 8y. 

The stability characteristics of the system with step input are 
also very similar to those of the system with zero input. 

The conclusion is that with y = 1.38 * 10 , the continuous-data 
system is always stable while the sampled-data system is stable for T 
less than 0.25 sec. 

For small sampling periods, the dynamic behavior of the sampled- 
data system is very similar to that of the continuous-data system. When 
T is large (but less than 0.25 sec) the overshoot of the step response 
of the sampled-data system becomes greater. However, the dynamic behavior 
of the sampled-data system may be improved by redesigning the controller. 

From this study it appears that a sampling period as high as 0.1 
second is feasible for the LST system. However, it should be noted that 
the conclusions are obtained with the existing system model. Other 
practical considerations such as noise, coupling effects and quantization 
errors, may restrict the sampling period to a lower value. 
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2. Design of the Continuous-Data LST System 


The purpose of this section is to carry out an optimal design of 
the continuous-data LST system. The strategy is that the continuous- 
data controller may be used as a basis for a digital redesign, and at 
the same time, a digital control can be designed using a completely 
independent approach. 

The system model of the LST was adopted~^n)m-refemice [1], and 
was later simplified from a 9th-order system to a 4th-order system in 
reference [2]. This simplification was justified from the standpoint of 
the system parameters, with no resulting loss of reality. 

The block diagram of the 4th-order LST system is shown in Figure 2-1. 
Since K q and K.| represent the parameters of the controller which are reported 
in reference [1], it is of interest to consider a complete redesign of 
the system. It was pointed out in [2] that with K Q = 5758.35 and = 
1371.02, the dominant CMG and vehicle modes are all with a damping ratio 
of approximately 0.707. However, in this report an attempt is made to 
arrive at a different control using the optimal control technique. 


2-1 . Decomposition of the LST System 

Figure 2-2 shows a state diagram of the system of Figure 2-1. It 
is important to note that the nonlinear loop of the CMG dynamics is valid 
only from a symbolic viewpoint. In other words, the diagram of Figure 2-2 
is obtained by treating N as a linear gain. For design purposes, the non- 
linear loop is deleted, and for computer simulation, the system diagram of 
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Figure JM’ should be used. 

Figure 2-2 also indicates that the control technique of state feed- 
back is used. In reality, the system of Figure 2-1 feeds back. two states 
in 0 y and 8 y only. 

The purpose of constructing the state diagram is so that we can represent 
the system in state variable form. The state equations of the system in 
Figure 2-2 with the nonlinear loop open are 

x(t) = Ax(t) + Bu(t) , (2-1) 

where 

0.1 0 0 

0 0 ~ 0 

J V 

0 0 0 V 

J G 
-K 

0 0 -K t -f- 

1 J G 

0 
0 
0 
K 

The control is given by 

u(t) = K 0 x(t) - K qXi ( t) - K^t) - K 2 x 3 (t) - K 3 x 4 (t) 







(2-4) 
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One of the advantages of using the state variable feedback is that 
the system can be designed in the sense of an optimal linear regulator. 
The performance index used for the optimization is 


J = 


[x. 1 (t)Qx(t) + u' (t)Ru(t)]dt 
o 


(2-5) 


where Q is a symmetric semi-positive definite matrix,, and R is symmetric 
and positive definite. The design objective is to determine the optimal 
control u(t) so that J in Equation (2-5) is a minimum, subject to the 
equality constraint of Equation (2-1). 

It is well known that the solution to this optimal control problem 


is 


u(t) = -R - 1 B ‘ Kx(t) 

where K is the solution of the algebraic Riccati equation. 
-KA - A'K + KBR'Vk - Q = 0 


( 2 - 6 ) 


(2-7); 


The solutions of the Riccati equation and the optimal control have 
been programmed on a digital computer. Table 2-1 gives the solutions of 
Kq, K ] , K 2 , and K^, and the corresponding eigenvalues of the closed-loop 
system when various weighting matrices Q are used, where 



( 2 - 8 ) 



Linear Regulator Design of the Continuous-Data LSI System (q 
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Since the states and are of primary interest, and are kept 
constant at 1.0 while q.j and q^ are varied. Also, R = 1. 

Several facts become clear from the results of Table 2-1. 

(1) With large values of q-j and that is, more weights on 
the states x-j and the feedback gains 1^ and become 
negligible. 

(2) The eigenvalue at -9700 is relatively insensitive to the 
various weighting matrices. 


2-3. Design by Eigenvalue Assignment and the Inverse Problem 


The development in the last section shows that it is difficult 
to have complete control of the eigenvalues of the closed-loop system 
by changing the elements of the weighting matrix Q. Since the original 
system from [1] with K q = 5758.35 and K-j = 1371.02 resulted in a rather 
good step response, it is interesting to find out if it corresponds to an 
optimal linear regulator solution. This question is known as the inverse 
regulator problem [3]. 

The state equation of Equation (2-1) should first by transformed into 
the phase-variable canonical form. Substituting the system parameters into 
Equation (2-2) yields 


0 1 



0 



6 x 10" 3 0 

0 0.4762 

-9700 -102.86 


(2-9) 
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O ' 

0 

( 2 - 10 ) 

0 ' 

9700 

The transformation which transforms Equation (2-1) into the phase- 
variable canonical form is 

10 0 O' - 

0 10 : 0 

x (2-11) 

0 0 0.006 0 

0 0 0 0.00286 

v ^27.742u " (2-12) 

The transformed state equation becomes 

i(t) = A^Jt) + B ]V (t) (2-13) 




where 



(2-14) 


(2-15) 
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The state feedback of the original system is described by 

u = -Gx = -[K q K } K 2 K 3 ]x (2-16) 

For the transformed system. 


0 
=n 

i — i 

1 

ii 

£ 

i 

it 

> 

H ] H 

2 H 3^ 


(2-17) 

Thus, using Equations 

(2-11), 

(2-16) and (2-17), G and H are related 

through 


’ 1 

0 

0 0 



c = H 

0 

1 

0 0 



u 27. 

,742 

0 

0 

0.006 0 





0 

0 

0 0.00286 




| 0.036 

0 

0 

0 



0 


0.036 

0 

0 


= H 

0 


0 

0.000216 

0 

(2-18) 


0 


0 

0 0.000103 j 


Thus, 







K Q = 0.036H ( 

3 





K ] = 0, 

.036H 1 








(2-19) 

K 2 = 0.000216H 2 





o 

II 

CO 

.000103H 3 





The inverse problem is 

that given the matrices A and B, and the 

feedback matrix G 

, is 

there a 

i positive definite R and nonnegative Q such 


that Equation (2-16) is the optimal control for the system of Equation 
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(2-1) with the performance index in Equation (2-5)? 

We form the characteristic equation of the closed-loop system in 
the 'phase-variable canonical form. 

p(s) = |sl - A ] + B-jH | = s 4 + (102.86 + H 3 )s 3 + (4623.7 + H 2 )s 2 

■ + H^s + H : (2-20) 

The characteristic equation of the open-loop system is 

ip(s) = | si - 7fc~| = s 4 + 102.86s 3 + 4623.7s 2 (2-21) 

It is shown in [3] that the feedback matrix H is indeed optimal 
for the choice of 

Q = D'D ' (2-22) 

and R = 1 , where 

' D = [d ] d 2 d 3 d 4 ] (2-23) 

The elements of D are the coefficients of the polynomial 
m(s) = d 4 s 3 + d 3 s 2 + d 2 s + d ] (2-24) 

where m(s) satisfies 

p(s)p(-s) = \p(s)ip(-s) + m(s)m(-s) (2-25) 

Substitution of Equations (2-20), (2-21), and (2-24) into Equation (2-25), 
we have the following relationships after simplification: 
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-1332.78 - = 2(H 2 + 4623.7) - (102.86 + H 3 ) 2 (2-26) 

(4623. 7) 2 + d 2 - 2d 2 d 4 = 2 H q - 2^(102.86 + H 3 ) + (H 2 + 4623. 7) 2 (2-27) 

2d ] d 3 - d 2 = 2 H q (H 2 + 4623.7) - H 2 (2-28) 

d 2 = H 2 (2-29) 

For 

K q = 5758.35, 

K 1 = 1371.02, 
k 2 = k 3 = 0, 


we have 


H q = 159954, 
H ] = 38083.9, 
h 2 = h 3 = 0. 


Thus, the last four equations lead to 

d ] = H q = 159954 (2-30) 

d 4 =0 (2-31) 

d 2 = 2H - 205. 72H, =-7.51 x 10 6 
3 0 1 


(2-32) 
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2 

Since is negative, we do not have a real solution for d^- Thus, the 
system with the prescribed feedback gains does not correspond to an 
optimal linear regulator solution. Equation (2-32) shows that to have 
a linear regulator solution, the following conditions must be satisfied 
(with l<2 = K3 = 0): 

H q > 102.86H-J (2-33) 

or 

K q > 1 02 . 86K-J (2-34) 

Eigenvalue Assignment 

An alternative to the linear state regulator method for designing 
linear feedback systems is the method of pole placement or eigenvalue 
assignment. In this method, the approach is to place the eigenvalues of 
the closed-loop system at certain desired locations by appropriate choice 
of the feedback gains. If the system is in the phase-variable canonical 
form, this method is directly applicable [3]. 

Consider the linear system represented by Equations (2-14), (2-15), 
and (2-17), the characteristic polynomial of this closed-loop system is 
as in Equation (2-20), 

| si - A ] + Bl H|= s 4 + (102.86 + H 3 )s 3 + (4623.7 + H 2 )s 2 

+ H-j s + H q (2-35) 

Let the desired location of the closed-loop eigenvalues be -a^ , -a 2 , 

-a^, -a^. The characteristic polynomial which yields these eigenvalues is 



. . J 

The desired feedback gains are obtained by equating the coefficients 
of the polynomials of Equations (2-35) and (2-36). Thus 

H o = a l a 2 a 3 a 4 

H 1 = a l a 2 a 3 + a l a 2 a 4 + a l a 3 a 4 + a 2 a 3 a 4 

^2 = a l a 2 + a l a 3 + a l a 4 + a 2 a 3 + a 2 a 4 + a 3 a 4 ~ 

H 3 = cxi + °i 2 + a .3 + <*4 - (2-37) 

In addition, it is desired that all elements of H be positive so that 
negative feedback is maintained. Once H is determined the feedback matrix 
G can be obtained from Equation (2-18). 

Table 2-2 shows the feedback gains G and H, obtained for several 


choices of eigenvalues. 
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3. Computer Simulation of the Simplified LST System with the Linear 
State Regulator and Eigenvalue Assignment Designs 

A computer simulation is presented here to show the response of 
the continuous-data LST system with the feedback gains designed by the 
linear state regulator method and the eigenvalue assignment method. 
Although the design was performed in Section 2 without the nonlinearity, 
the simulation is of the nonlinear system. As mentioned earlier, the 
system of Figure 2-1 is used for this purpose. 

The numerical data for the system and nonlinearity are 

J v = 10 5 

J G = 2J . 

K p = 216 

Kj = 9700 
H = 600 

y = 1 .38 x 10 5 



For the computer simulation, the input to the LST system x is set to 

O O , 

5 x 10 K so that the final value of 0,, = 5 x 10 . ATT initial conditions; 

o V 

are zero. The following quantities are plotted 

0y = vehicle position (radians) 

. ojy = vehicle velocity (radians/second) 

0 G = Gimbal position (radians) 



uj_ = Gimbal velocity (radians/second) 
b 

TgP = Torque output of the nonlinearity (ft. lb.) 

Error = Error input to the CMG. 

= X - K q 0 v - K-jlOy 

In all the simulations, only the states 0y and are fed back: 
Although the designs yielded feedback from all four states, the gains 
« 2 , K 3 are ignored as they are very small and are also not accessible 
in the system which is simulated. 

Figures 3-1 and 3-2 show the results with a design using the eigen^ 
value assignment method. The feedback gains are K q = 5773, K, = 

2032. The corresponding eigenvalue locations are -500, -40, -4 ± j4. 

Figures 3-3 and 3-4 show the results with a state regulator design. 
The feedback gains are K Q = 10000, =6579.5. The system eigenvalues 

are -9700, -3.08, -1.545 ± j2 . 626 . The weighting matrices used for this .. 

o 4 

design are Q = diagonal [10 10 1 1], R = 1. 

Figures 3-5 and 3-6 show the results with another state regulator 
design. Here « 0 = 7071 , K-j = 5220, the eigenvalues are -9700, -2.75, 

-1 .375 ± j 2 . 33 and the weighting matrices are Q = diagonal [5 x 10 7 5000 

1 1], R = 1. 

The results of Figures 3-1 through 3-6 show that the linear design 
methods of Section 2 can yield closed-loop systems with very acceptable 
response characteristics. In comparison with, the system. used in reference 
[2], the new designs result with no overshoot in the response of 0y. The 
state regulator designs are, however, slower in comparison to the eigen- 
value assignment design and the system in reference [2]. 



e G (*io' J ) a> v (*io‘ 6 ) e v (xio'°) 

0.00 ’ 0.01 0.02 0.00 0.08 0.16 0.00 0.025 0.05 


34 


System 1, y = 1. 38*1 0 5 , No Sampler 
Pole Placement Design 




TIME 



TIME 


Figure 3-1 





1 ' 

0.00 

1 1 
1.50 

n 1 r 

d.OO 4.50 

1 1 ME 

r 

6.00 

n ^ 
7.50 

1 

9.00 



Fiqure 3-3 





TiI'O- 




0 G (xlO _H ) ^(xlO -0 ) 0 v (xlO"°) 

U.00 ’ U.D8 0.16 U.00 U.04 0.08 0.00 0.025 0.05 


38 


System T, Y = 1 .38*10 5 , No Sampler 
State Regulator Design 



TIME 



1 IME 



TIME 


Figure 3-5 




Error(xlO"^) ^SF^® ) ojg( x 10 ) 

- 0.05 0.15 0.35 -< 4.23 5.77 15.77 - 0.06 0.09 0.24 


System 1, Y = 1.38*10^, No Sampler 
State Regulator Design 



0.00 


T 

1.50 


“1 

3.00 


T 


TIME 


n — 

4.50 


1 1 

6.00 7.50 


9.00 



TIME 



Figure 3-6 


40 


4. Digital Redesign of the Large Space Telescope (LST) System 

In this chapter, the simplified digital LST system is designed 
by use of the point-by-point method of digital redesign [6], Three 
different continuous-data control laws are considered. 

Control Law A 

This is the original control law of the LST system, obtained by 
classical techniques for a damping ratio of 0.707. The gain matrices 
for this case are 

G(Q) = [5758 1371 0 0] 

E(0) = [5758] 

This control has been used extensively in analyzing the stability 
of the LST system with the nonlinear CMG gimbal friction [2]. 

Control Law B 

This control law is obtained by the eigenvalue assignment method [7]. 

G( 0) = [5773 2032 0.97 0.04] 

E(0) = [5773] ' . . 

Control Law C 

V 

This control law is obtained by the linear regulator design method [7]. 
G(0) = [7071 5220 10.6 0.99] 


E(0) = [7071] 
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The simplified model of the continuous-data LST system is used 
in its decomposed form [7]. 

xU) = Ax^(t) + Bu(t) (4-1) 

where 

0 10 0 

0 0^-0 

v 

0 0 0 ~ 

J G 

-K 

0 0 -k t -£ 

1 J G 

0 
0 
0 
K 

The control is given by (with zero input) 

u (t) = -G(0)x(t) 

The numerical values used are 

H = 600 
J v ■ 1 x lo 5 

J G =2.1 

K p = 216 

Kj = 9700 
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It should be noted that the system model of Equations (4-1), 

(4-2). and (4-3) is obtained by neglecting the nonlinearity of the CMG. 

This is necessary in order to obtain a state variable representation 
of the system. In reality, however, only the first two states of the 
system, 0y and 0^, will be fed back. Thus, the last two gains of control 
laws B and C are for redesign purposes only; they will not be used in 
the simulations. 

The closed-loop eigenvalues of the LST system, with the three 
continuous-data control laws are: 

Control law A: s^ = -46.86 + j39.04 

s 2 = -46.86 - J39.04 
s 3 = - 4.57 + j4. 70 
s 4 = - 4.57 - j4.70 

Control law B: s^ = -471 

s 2 = -11.1 
s 3 = -3.97 + j'3.85 
= -3.97 - j 3 . 85 

Control law C: s^ = -9700 

s 2 = -2.77 

s. = -1.37 + j 2 .32 

s 4 = -1.37 - j2.32 1 

The point-by-point method of state matching is used to digitally 
redesign the LST system for each of the continuous -data control laws. 
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With each control law, a weighting matrix H is chosen and the redesign 
is performed for a range of sampling periods T. With the feedback gain 
matrix G W (T), the digital system is represented by 


x(k + 1) =■ <MT)x(k) •+ 8(T)u(k) 

(4-4) 

u(k) = -G w (T)x(k) 

(4-5) 

where 


4>(T) = e AT 

(4-6) 

f T At 

6 (T) = e HT dtB 

(4-7) 

'o 



The z-plane characteristic roots of the sampled-data system are given 
by the eigenvalues of (4>(T) - e(T)G w (T)). These roots are calculated for 
each value of T with some very interesting results. 

The results with control law A are shown in Table 4-1 for H = 

[1 1 0 0], Table 4-2 for H = [1 1 1 1], and Table 4-3 for 

H = [1 0 0 0]. The results with control law B for H = [1 1 

0 0] are shown in Table 4-4, and the results with control law C 

for H = [1 1 0 0] and H = [1 0 0 0] are shown in Tables 

4-5 and 4-6, respectively. 

The z-plane characteristic roots of the redesigned systems in 
Tables 4-1 through 4-6 are shown in Figures 4-1 through 4-6, respectively. 

These results show that a large variety of sampled-data systems 
are available to control the digital LST. It is apparent that the 
choice of H plays a very dominant role in the resulting gain G^(T). 

In fact, with control law A, while one choice of H (Table 4-1 and 
Figure 4-1) provides a stable sampled-data system for a wide range 
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of sampling periods T, another choice of H (Table 4-2 and Figure 4-2) 
yields a system which is unstable at higher sampling periods (T greater 
than 0.5 sec). Still a third choice of H (Table 4-3 and Figure 4-3) 
yields a system which is stable only at very low sampling periods 
(T = 0.06 sec or less). With control law C, and H = [1 1 0 0] 

or H = [1 0 0 0], the redesigned system is unstable for very small 

sampling periods (T less than 0.018 sec and 0.02 sec, respectively). 

In general, it appears that control laws A and B are more effective 
when digitally redesigned. 

Some computer simulations are now presented to show the effects 
of digital redesign with the various control laws. As before, only the 
first two feedback gains are utilized, and the simplified LST system is 
simulated on a digital computer with the CMG nonlinearity included. The 
parameters of the nonlinearity are y = 1.38 x lO'*, Tgpg = 0.1. In all 

O 

the simulations, the input x = 0 and the initial value of 0 V (O) = 1 * 10 , 

with all other initial states equal to zero. 

The following simulations have been performed: 


Figure No. 

Control 

Law 

H 

Sampling Period 
T (sec) 

4-7, 4-8 

A 


- 


Continuous-data System 

4-9, 4-10 

A 

[1 

1 0 

0] 

0.02 

4-11, 4-12 

A 

[1 

1 0 

0] 

0.1 

4-13, 4-14 

B 


- 


Continuous-data System. 

4-15, 4-16 

B 

[1 

1 0 

0] 

0.02 

4-17, 4-18 

B 

[1 

1 0 

0] 

0.1 
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Figure No. 

Control 

Law 

H 

Sampling Period 
T (sec) 

4-19, 4-20 

C 

- 

Continuous-data System 

4-21, 2-22 

C 

[1 10 0] 

0.02 

4-23, 4-24 

C 

[110 0] 

0.1 


In each simulation, the following quantities are plotted: 

e v = vehicle position (radians) 

w v = vehicle velocity (radians/sec) 

0g = Gimbal position (radians) 

(d g = Gimbal velocity (radians/sec) 

Tgp = Nonlinearity Torque (ft-lb) 

Error = x - K q 0 v - 

= Error input to CMG controller. 

The simulation results show that adequate digital control schemes 
can be obtained for wide ranges of sampling periods by appropriate re- 
design of the feedback gains. Again, it appears that the method of 
redesign is more effective in the case of control laws A and B. 
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Y = 1.38, T = 0.1 sec 
Control Law A 



ERROR (x 







ERROR ( X 1 0 °) t GF^ x1 ° > w G (*10 *) 

-0.02 -0.01 0.01 j -1.04 -1.51 0.96 -0.02 -0.002 0.016 
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Y = 1.38, T = 0.02 sec Control Law B 




0.00 1.50 - -3.00 4.50 , - 6.00 ... 7.50 9.00 

TIME 



Figure 4-16 



- 0.02 





n ' i 1 1 ■ 1 1 

3.Q0 4.50 6.00 7.50 9.00 

TIME 


1 1 i 1 1 1 1 1 

10 4.50 6.00 7.50 9.00 

TIME 


Figure 4-20 










' I 

0.00 

1 I 

1.50 

i ' l T 

3.00 4.50 

TIME 

I ' 

6.00 

l ; ; r 
...7.50 

9.00 



Figure 4-23 
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5. Stability Considerations and Constraints on the Selection of the 
Weighting Matrix of the Digital Redesign Technique 

It has been reported [6] that given a continuous-data system 

x(t) = [A - BG(0)]x(t) (5-1) 

*. t ; ' 

the solution of the feedback matrix G(t) of an equivalent digital 
system, designed with the point-by-point state comparison method, 
must satisfy the following equation: 

0(T)G(T) ■ e AT - e CA-BGC0)lT (5 . 2 ) 

Since 0(T) is usually not square, we cannot solve for G(T) directly from 
the last equation. One remedy to the problem is to introduce a weighting 
matrix H, such that the inverse of H0(T) exists. Then, 

G W (T) = [H0(T)] _1 H[e AT - e [A ’ BG(0)]T ] (5-3) 

However, the weighting matrix H cannot be chosen arbitrarily. The 
solution in Eq. (5-3) is significant only if the digitally redesigned : 
system is stable. 

In chapter 4, it has been demonstrated in the digital redesign, 
of the LST system that for some sampling period T and some H, the 
resultant G W (T ) gives rise to an unstable closed-loop digital system. 

This means that given the continuous-data control system, the weighting 
matrix H cannot by chosen arbitrarily. The conclusion is that if the 
closed-loop digital system is unstable, the solution to G W (T), corresponding 
to the selected H, will be meaningless. 
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The problem now is to find the condition under which an H can 
be found such that the digital system is stable. 

Stability of the Closed-Loop Digital System 

The state equations of the digital system are written 


x[(k + 1)T] = (j)(T)x(kT) + 6(T)u(kT) 

(5-4) 

where 



<t>(T) = e AT 

(5-5) 

e(T) = 

f T e AX dAB 

(5-6) 


’o 


The feedback control is 
u(kT) = -G(T)x(kT) 

(5-7) 

Then, En- (5-4) becomes 


x[(k + DT] = [<h(T ) - 0(T)G(T)]x(kT) 

(5-8) 


The digital system is stable if all the eigenvalues of [4>(T) - 0(T)G(T)] 
are located inside the unit circle |z| = 1. Since $(T) and 0 (T ) are 
known once the sampling period T is specified, the conditions on G(T) 
for stability can be established using well-established techniques. 

Let - 

e At _ e [A-BG(0)]T ; D(T) " (5 

and premultiplying both sides of Eq. (5-2) by the 1 x n matrix H, 
we have r;c . - „ 



/ 
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H8(T)G w (T) = HD (T ) 


(5-10) 


where G (T ) has been replaced by G W (T) to indicate the Weighed matching 
of states. 

Taking the matrix trafts pose ? ’6nboth si' des' , of , the - Tast^equation , 
we get, . , ' • 


G^(T)e'(T)H' = D'(T)H‘ 


(5-11) 


Rearranging, Eq. (5-11) becomes 


[Gy(T)9'(T) - D' (T)]H' = 0 


(5-12) 


This equation represents a set of n linear homo'geneous equations which 
have nontrival solutions if and only if the following condition is 
satisfied: 


Gy (T ) 0 1 (T) - D f (T)| = 0 


(5-13) 


which is also equivalent to 


|e(T)G w (T) - D(T) | = 0 


(5-14) 


Thus, if Eq. (5-14) is satisfied, there always a nonzero H which 
will satisfy / 


G W (T) = [ H 0 ( T ) ],~ 1 HD (T y 


(5-15) 


Illustrative Example 


Consider the continuous-data system 


• '/ ' 


I 
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x(t) : 'Ax(t) t Bu(t) 
u(t) = -G(0)x(t) 

• ... ,‘\ _ . t 

where 



' 0 

1 ' 


O’ 

A = 


■ 

B' - . 



0 

0 


1 


G(0) = [2 3] 


(5-16) 

(5-17) 



It is desired to design a digital system which will match the response 
of the continuous-data system at the sampling instants. The sampling 
period is 1 second. 

The following. matrices are computed: t 


<fr(T) * 


• 1 

T ' 


' 1 

1 

0 

1 


0 

1 . 


0(T) = 


0.5 ' 
1 


D(T) 


0.767 1.233 

0.465 1.097 


The characteristic equation of the closed-loop system is 


F(z) = | zl - <(>(T) + 6G W (T) | = z 2 + [-2 + O.SG^T) + G 2 (T)]z 

+ [1 - 0. 5G ] (T) - G 2 (T) + TG 1 (T)] = 0 

(5-18) 


where G^(T) and G 2 (T) are the elements of G^(T); that is. 
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G W (T) = [G 1 (T) G 2 (T)j (5-19) 

Using the Schur-Cohen stability criterion/ the roots of Eqi (5-18) 
are all inside the unit circle if 

F ( 0 ) = 1 + 0. 5G -j (T ) - G 2 (T) <1 .... .. (5-20) 

F ( 1 ) = G 1 (T) >0 (5-21) 

F ( - 1 ) = 4 - 2G 2 (T) > 0 (5-22) 

These conditions on G^(T) and G 2 (T) are plotted in the parameter 
plane of G 2 (T ) versus G^(T), as shown in Fig; 5-1. •••; 

Having established the conditions. on the elements of G^(T) for the 
stability of the digital system, we turn to the condition under which 
an H exists which also satisfies Eq. (5-15). 

Equation (5-14) leads to 

0.5G-J (T) - 0.767 0.5G 2 (T) - 1.233 

|0(T)G W (T) - D(T)| = =0 

w G 1 (T) - 0.465 G 2 (T) - 1.097 

(5-23) 

or 

• -0,535G 2 (T) + 0^684G.|(T) + 0.268 = 0 (5-24) 

Equation (5-24) represents a straight line in the G 2 (T) versus 
G^T) parameter plane. The intersect between the line represent 
by Eq. (.5-24) and the stable region gives the stable trajectory for 
G}(T) and G 2 (T), as shown in Fig. 5-1. 



81 


If the intersect between Eq. (5-24) and the stability region of 
G-j (T) and (^(T) is convex jn general , the vertices of the intersect 
can be used to find the, bounds .on the weighting matrix H. 

In the present case, the vertices of G-j (T) and G,,(T) are at 
(0, 0.5023) and (1.171, 2). 

Substituting the vertices of .G-j(T) and G 2 (T) in Eq. (5-12), 
we have the two boundary equations for the elements of H = [h^ h^]. 


G-j (T ) = 0, 

h| = -0.606h 2 

(5-25) 

G ] (T) = 1.171 

h ] = 3.875h 2 

(5-26) 


Figure 5-2 shows the region in which h^ and h ? should lie so that 
Eq.' (5-3) will always yield a state .feedback control such that, the 
digital system is stable; ''' - r 






Figure 5rl . Stable region and locus for nontrival solution of H. 





Figure 5-2. Allowable regions of-h, and hi for stable digital system. 
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6. Real izatioW of ^State Feedback by • 6ynam:i£’ /Co^ray^’?.'--.? 

One of 'the unique characteristics of modern control theory is that 
optimal control is often realized by state feedback. For instance, it is 
well known that if a system is completely controllable, its eigenvalues 
can be arbitrarily assigned through state feedback, and the optimal 
linear regulator design always leads to a state feedback solution. Un- 
fortunately, in practice, not all the state variables of a physical 
system are accessible. Considerable amount of results have been reported 
in the past on the design of optimal systems with partial state feedback. 

The basis of the classical control system design^is' that the 
configuration of the controller is selected a priori. The controller 
used in practical systems usually assume the form of cascade or feed- 
back controllers, or a combination of these. In these cases, only the 
outputs of the system are fed back. One advantage of the classical 
controllers is that they can be implemented often by passive filters or 
electronic circuits. . 

In this, chapter, we shall present a method whereby a system with 
state feedback is approximated by a system with a cascade controller. 

How the state feedback is determined is immaterial for the present 
analysis; it! can be obtained from the pole-location solution or the 
Riccati equation solution, or some other optimal control design methods. 

Continuous-Data Systems 


Consider the system 
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x(t) = Ax(t) + Bu(t) 


( 6 - 1 ) 


J'('t) = Cx(t) + Du(t) 


( 6 - 2 ) 


where 


x_(t) = n x 1; state vector 
u(t) = r x 1 input vector 
y(t) - in x 1 output vector 




A, B, C, and D are coefficient matrices of appropriate dimensions. 
Assume that state feedback is given such that 


u(t) 'K-Gx(t)’ 


(6-3) 


where G is an r x n feedback gain matrix. 

\ The design objective is to approximate the system of Fig. 6-la 

which is described by Eqs. (6-1), (6-2), and (6-3), by the system of 

\ 

Fig. 6-lb which has a feedback controller with feedback from the out- 
put Variables. Let the transfer relation of the controller be repre- 
sented x by 

U(s)\-H(s)Y(s) 


(6-4) 


where H(s) is f^e controller transfer function matrix: 


H(s) = 


f H 11 (s x ). H 12 (s) 

H 91 (s) H„(s) 


'21 


22 ' 


r2 


H 1m< s > 

V s ) 


H r] (s) H„,(s) - . J H_(s) 


rm' 


(6-5) 
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Let H . . ( s ) (i = 1, 2, r; j = 1, 2, .... m) be a pth order 

* *J 


transfer function. 






K . . (-1- + a'. '.-is '+- a.. . • '•-+ a. - -S P ) . 

t-i' "ill it/ iin • 


jjr ij2 


2JP 


( ] + 6 ijl S+ ®ij2 S + -•* + e i3P sP > 


( 6 - 6 ) 


The transfer function H..(s) is expanded into a Taylor series 

* J 


about s = 0, 


H U (S) ■ K fo J 0 d ijk ! 


... 18:7,) 


where 


d ijk 


3 V s > 


8S 


s =. 0 


(6-8) 


Evaluating the coefficients of .the Taylor series, we have 


d ijo - ’ 


d 1jl ° “ijl ' B ijl 


and for k > 1 , 


k-1 

d ijk " a ijk ” e ijk " y ^ 6 i j(k-v) d i jv 


(6-9) 


Since the state feedback represents the feedback of the system 
output and its higher-order derivatives, a truncated series expansion 
H..(s) may be used as a dynamic implementation of state feedback by 

* J 
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feeding back only the output variables. 

If the infinite series of (6-7) converges, we may approximate it 
truncating it after p terms,, where p is not. yet specified. Let us 
introduce the. notation, H.. (s) , for the truncated version of H..(s); 

1 J P * J 

then 


V s! "Si X d iik sk 


( 6 - 10 ) 


where d.j k is as defined in Eqs. (6-9). 

Substituting Eq. (6-10) in Eq. (6-4) for the elements of H(s), 
we have 


“ill <Wl)3 ^2^120 d 12! ••• d 12(p-1)l ••• d lml - Vp~l> ] j 

V,( S ) 

' i 

^,(s) 

^1^210 d 211 ••• d 21(p-1) J *22 Cd 220 d 221 - d 22(p-l)l •••' S, [d 2n,0 d 2n, - d 2n(p-l,l 

!• 

n 1 


s p *\(s) 

- * - - 'j 

» 

; • | 



sY 2 (s) 


. ! 

s p '’’v 2 's) ; 

f 

i 

i 

i 

j.: j 

J 

j V‘> ! 

i 

i 

i 

Sir's) • 

m • 

> i 

K r1 Ld rlO. d HT ••• d rl(p-l)l l ‘rtf d r20 d r 21 — d r2(p-l)l — K J*roO d nnl ”• d rm(p-t) J j [ S ”'V 51 j 


( 6 - 11 ) 
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The elements of the last equation are rearranged to give. 
t' K 11 d no *12 d 120 ••• K 12 d 121 ••• t*11 d n(p-l) *12 d 12(p-1) 

C *21 d 210 K 22 d 220 "• lt 2m d 2mO^ < 21 d 2n *22 d 221 K 2m d 2ral^ J " flC 21 d 21(p-V) *22 d 22(p-l) 


K lro d lin(p-1 p 
*2m d 2»>(p-l : 


Xts) 


*V(s> l 

, 

s ? y(s! !j 


i 


^r! d rH) X r 2 d r 26 •" K rT! d rm0^ l< rI d r11 lt r2 d r21 *• 


K n» d r«n 3 •” ^rl d rl(p-1) } K r2 d r2(p-1) *” Vi*!p-1 )l f 


s p-V XCs) 


( 6 - 12 ) 


i.ii 


The time-domain equivalence' of the last equation is 


f z(t) 


i(t) 


u(t) 



lit) 


(6-13) 


. x <p_ 1 ) (t> . 

where F denotes the r * mp coefficient matrix in Eq. ( 6 - 12 ). 

From Eq. ( 6 - 2 ), 

/(t) = C_x ( t ) + Du(t) 

= (C - DG)x(t) ( 6 - 14 ) 


Then 
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i(t) = (C - DG)(A - BG)x(t) 


= (C - DG) (A - BG) P_1 x(t) 


Substituting the last three equations in Eq. (6-13), we have 


u ( t ) = -F 


(G t DG) 

(C - DG) (A - BG) 


x(t) 


(C - DG) (A - BG) P_1 


Comparing Eq. (6-17) with Eq. (6-2), we have 


(C - DG) 

(C - DG) (A - BG) 



l (C - DG)(A - BG ) P_1 j 
(r ,x jnp) (mp x n) (r x n ) 


In order to solve for F from the last equation, we write 


(6-15) 


(6-16) 


(6-17) 


(6-18) 
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(C - D6) 

(C - DG)(A - BG) 

(C - DG) (A - BG) P_1 

if mp = ri, or p = n/m. This means that if n/m is an integer, we may 
truncate the Taylor series expansions of G..(s), i = U 2, r, 

* J 

j = 1, 2, .... m, at p = n/m terms. If n/m is not ah integer, we may 
choose p to be an integer which satisfies 

- < p < - + 1 (6-20) 

m K m 

Since F is r x mp, there will be (r)(m)(p) unknowns. However, there 
are only rn equations in Eq. (6-18). Thus, r(mp - n) of the elements 
of F may be assigned arbitrarily. 

The solution of F from Eq. (6-19) also depends on the existence of 
the inverse in the equation. 

It should be noted that solution of the elements of F gives only 
the values of the coefficients in Eq. (6-7). The coefficients of the 
transfer function of (6-6) still have to be determined using Eq. (6-9)'. 

In general, there are more unknowns than equations in Eq. (6-9). This 
simply means that in the ideal situation we simply set 



i d ./. . = oU , • » ‘ • 

i j k ijk 

and all 6.. = 0, for k = 1, 2, ..., m. However, for a physically 
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realizable transfer function, H.-. (s) must not have mor?e zeros than 

* Jr 

poles. Therefore, the values of should be assigned such that 

the dynamic behavior of the overall system is not appreciably affected 

by the presence of 3^^, k = 1 , 2, . . . , m. This is similar to the 

classical design practice of designing the zeros of H.. (s) to control 

i ' 1 JP 

the dynamic behavior of the system, while placing the poles of H^ p (s) 
so that, they do_not have appreciable effects on the system performance. 

Single-Variable Continuous-Data Systems 

When the control u(t) is a scalar, u(t) = -Gx(t), where 


G =, [g 1 g 2 ... g n ] 


N 


Eq. (6-6) becomes 


( 6 - 21 ) 


H(s) = 


K(1 + cx-jS + o^s + 


(1 + B-|S + B 2 s^ + 


v"> 


t 8 n S "> 


( 6 - 22 ) 


Then, 


H(s) = K ( 1 + d^ + d 2 s + 




(6-23) 


where 


^lc-1 

d k * “k.--' 6 k - j, 6 k-v d » 


(6-24) 


for k = 1 , 2, . . . , n-1 . 

Equation (6-19) becomes 
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F = K[1 d 1 d 2 


' C - DG 

(C - DG) (A - BG) 


<W = G 


(C - DG) (A - BG) 0 ' 1 


•(6-25) 


Single-Variable Continuous-Data Systems in Phase-Variable Canonical Form 

If the system to be controlled is in the phase-variable canonical 
form, then 


A = 


0 1 0 ... 0 

0 0 1 ... 0 


0 0 0 ... 1 


-a. 


-a 


n-1 


~ a n-2 



% 

(6-26) 


B = 


0 ' 
0 


(6-27) 


and the output equation if characterized by D = 0, and 

C = [1 0 0 ... 0] (6-28) 

the formulation given in the preceding section is further simplified. 

Since £ = 0, and CBG =0, 
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■ C - DG 


c 

(C - DG) (A - BG) 

_ 

CA 

(C - OG) (A - BG)" -1 


. ca"' 1 . 


= I (identity matrix) 


(6-29) 


Then, Eq. (6-25) becomes 


F = K[1 d 1 d 2 ... d n-1 ] = G 


■ [g l g 2 




(6-30) 


and 


K = g. 


“1 = 


_ 9 1 3 1 + 9 2 


(6-31) 


a k = 


d k+l ^ 

g| + e k + v ^ e k-v 9 v+l 


k = 1, 2, ... . , n. 

Equivalent Cascade Controller 

The development carried out in the preceding sections is based on 
a controller being placed in the feedback path of the system as shown 
in Fig. 5-lb. When the reference input r(t) is zero, that is, when the 
system is a regulator, it does not matter whether the controller H(s) 
is in the forward path or the feedback path. However, when the input 
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is not zero, it may be desirable to determine an equivalent system which 
has the controller in the forward path of the system as shown in Fig. 5-2. 

In the following, single-variable notation is used for simplicity. The 
problem is t'b find the transfer function of the cascade controller G c (s) 
so that the closed-loop transfer functions of the two systems with feed- 
back controller and the forward-path controller are identical. The 
solution of G c (s) is 

G C (s) = T V G(s)[H(V) ~T (6-32) 

where 

G(s) = C(sl - A)"]b D = 0 (6-33) 

In general, given G(s), and having determined H(s), the order of 
G c ( s ) will usually be higher than that of H(s). 

The following example will illustrate the design method outlined 
in the preceding sections. 


Example 6-1 

Consider that the dynamic equations of a linear time-invariant 
system are given by 


x(t) = Ax.(t) + Bu(t) 
y(t) = Cx(t) 


(6-34) 


■ 0 

1 

0 


' 0 ' 

0 

0 

1 

B = 

0 

0 

0 

-1 


1 


where 
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C .« [1 0 0] 

Since the system is completely controllable, we may assign the 
eigenvalues of the system arbitrarily. Furthermore, the state equations 
are already in phase-variable canonical form. With state feedback, 
u = -Gx, the closed-loop transfer function of the system is 



1 

s 3 + (g 3 + 1 )s 2 + g 2 s + g ] 


(6-35) 


The characteristic equation is 

s 3 + (g 3 + Ds 2 + g 2 s + g] = 0 (6-36) 

Let us assume that we wish to place the eigenvalues of the closed- 
loop system at s = -10,-1 + jl, and -1 - jl . Then, Eq. (5-36) gives 

9-| - 20, g 2 - 22, g 3 - 11 


or 

G = [20 22 11] (6-37) 

Now consider that the states x 2 and x 3 are not directly accessible, 
and it is desired to approximate the state-feedback solution by a 
feedback controller and output feedback. Since the system is of the 
third order, n = 3, the dynamic controller may be of the second order; 
that is, 

2 

1 + a, s + a„s 

H(s) = K ^ (6-38) 

1 + 3-jS + B 2 s 
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In the present case, the results of Eq. (6-30) may be used. 


Thu Si 


K = g ] = 20 
g i e i + g 2 'i 

a, ■ -L-L- 1 ' * 6, 


a 2 = ■£+ h + e i ^ = °- 55 + 1 - 1e i + S 2 

I I; / 


Assuming that physical circuit elements allow the selection of 
B-j and B 2 to be relatively small as compared with the resulting values 
of a-| and o^, we let 6-j =0.15 and 6^ = 0-005. Then, 


= 1.25 
02 = 0.72 


The transfer function of the feedback controller is 


H(s) = 2880 ■?— t .lr.. 7 ,36s. ±.l.-.3.89 
s + 30s + 200 


The closed-loop transfer function of the system with the cascade 
controller is 


Yf si _ _ 2880(s 2 + 1.736s + 1.389) 

si 5 + 31s 4 + 230s 3 + 3080s 2 + 5000s + 


4000 


A comparison of the step responses of the system with state 
feedback and the system with the feedback controller is shown in 
Fig. 6-3. 


(6-39) 


(6-40) 




Figure 6-1a. Block diagram of system with state feedback. 




Figure 6-lb. Block diagram of system with cascade controller 
from output feedback. . 
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Discrete-Data Systems 

The dynamic controller design technique described in the last 
section can be applied to discrete-data systems. Consider the dynamic 
equations, 

x(k +1) =. Ax(k) + Bu(k) 

X(k) = Cx(k) + Du(k) 
where 

x(k) = n x 1 state vector 
ujk) = r x 1 input vector 
yjk) = m x 1 output vector 

A, B, C, and D are coefficient matrices of appropriate dimensions. 
Assume that the state feedback is used such that 

u(k) = -Gx(k) 

where G is an r * n feedback gain matrix. 

Let the controller be modeled as a feedback controller with the 
transfer function relation, 

U(z) = -H(z)Y(z) 

where H(z) is given by 


•H„(z) 

H ]2 (z) 

• • • H lm^ z ^ 

H 21 (z) 

H 22 (z) 

. . . H 2m ( z ) 

H r] (z) 

H (Z) 
rr' ' 

• • • H ™< z > 


(6-41) ‘ 
(6-42) 


(6-43) 


(6-44) 


H(z) = 


(6-45) 
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Let H.-(z) be a pth-order transfer function, 

I J 

K-.fz ' 3 + a..,z^ ^ + a. ^ + + a.. ) 

H ( Z ) = _uJ iiJ iii_ lie! 

i j ' ,P . R 7 P _1 + « 7 P -2 + +« 

z + 6 ijl Z + 0 ij2 Z + ••• + & ijp 


•Let us expand H.,(z) into a Laurent's series about z = 0, 

I J 


V 2 > ■ K ij J 0 d ijk z ‘ 


where 


(6-46) 


(6-47) 



Similar to the development in Eqs. (6-1 1 ) and (6-12), the time- 
domain correspondence of Eq. (6-44) is 

zOO 

Z(k - 1) 

u(k) = -F (6-50) 

. Z(k - P + 1) . 
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where F is identical to the r * mp coefficient matrix defined in 
Eqs. (6-12) and (6-13), except with its elements correspond to the 
coefficients of Eq. (6-46). 

From Eqs. (6-42) and (6-43), 

l{k) = (C - DG)x(k) (6-51) 

Thus, 

£(k - 1) = (C - DG)x(k - 1) (6-52) 

Also, 

Ax(k - 1) = x(k) - Bu(k - 1) 


= x(k) - BGx(k - 1) (6-53) 

Therefore, 

x(k - 1) = (A - BG) _1 x(k) (6-54.) 

Substitution of Eq. (6-54) in Eq. (6-52), we have 

X(k - 1) = (C - DG)(A - BG) -1 x(k) (6-55) 

Similarly, 

X(k - 2) = (C - DG) (A - BG)' 2 x(k) (6-56) 

lik - p + 1) = (C - DG) (A - BG)' p+1 x(k) (6-57) 

Thus, Eq. (6-50) becomes 
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U(k) 


' C - DG 

(C - DG) (A - BG )" 1 


l (C - DG) (A - BG ) _P+1 


x( tO- 


Companng Eq. (6^58) with Eq. (6-42), we have 


' C - DG 

(C - DG)(A - BG )' 1 
(C - DG) (A - BG )' P+1 


if mp = n, or p = n/m, and the indicated inverse exists. 

For a single-input, single-output system, r = m = 1. 
if D = 0, Eq. (6-59) is simplified to 


F = K[d o d ] d 2 ... d nl ] = G 


Example 6-2 

Consider the sampled-data process shown in Fig. 6-4. 


CA 

CA* 


-1 


CA 


-n+1 


(6-58) 


(6-59) 


Furthermore, 


(6-60) 


The z-transfer 


function of the process is 
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G(z)-M- °-3W^ 0.264 

g ' 2) z - 1.368Z + 0.368 

The state equations of the system can be written in the form of Eq. 
with 


(6-61) 

(6-41) 



0 

1 ' 


■ o • 

A = 



B = 



. -0.368 

1.368 . 


1 


The output equation, Eq. (6-42), is 

y(k) = Cx(k) = [0.264 0.368]x(k) (6-62) 

Let the state feedback be denoted by 

ii(k) = -Gx(k) = -[g 1 g 2 ]x(k) (6-63) 

The characteristic equation of the closed-loop system is written 
| (zl - A + BG) | = z 2 + (g 2 - 1.368)z + (g ] + 0.368) 

Let us select the feedback gains as g^ = 0.132 and g 2 = 0.368 so that 
the eigenvalues are at 

A.| = 0.5 + j0.5, = 0-5 - j0.5 

To obtain an equivalent cascade controller to replace the state 
feedback, we let 

H(z) = K = K (1 (6-64) 
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Thus, Eq. (6-60) gives 


F = K[1 d^ = G 


f C 
CA 


-1 


-1 


[0.132 0.368] 


0.264 

1.35 


= 0 . 86[ 1 -0.0708] 


0.368 

-0.717 


,-1 


(6-65) 


Then, K = 0.86, d ] = -0.0708. From Eq. (6-48), 


d i ■ “i - B i 


Selecting 6-j = 0.0005, we have = -0.0703. The transfer function of 
the feedback controller is 


H(z) = 0,86 


z - 0.0703 
z - 0.0005 


( 6 - 66 ) 


The overall system is shown in Fig. 6-5. 



(t) U*(t) 

z.o.h. 


1 

T * 

T = 1 sec 


s(s + 1) 


G(z) 


Figure 6-4. A sampled-data process. 



0.86 


Z-0.07CH 


z-0.000a 


y(k) 




Figure 6-5. Closed-loop sampled-data system. 
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7. A Numerical Technique for Predicting Self-Sustained Oscillations 

i 

in the Nonlinear LST System With the Continuous and .Discrete 
Describing Function Methods \ . . . 


7.1 Introduction 


\ 


It has been demonstrated [2] that the continuous and discrete 

i 1 ' ’ 

describing function methods are useful tools for predicting the 
existence of self-sustained oscillations due to the nonlinear CMG 

i w \ 

friction in a single-axis model of (the LST. The corresponding stability 
equations which have to be solved are 


1 + N(A)G(jw) = 0 


(continuous case) 


and 


1 + N(A, n)G(T, n) = 0 


(discrete case) 


(7-1) 


(7-2) 


In Reference [2], these equations were solved graphically, and the 

i 

conditions for existence of selfr sustained oscillations were established. 

/ 

/ 

In each case the intersection of the -1/N curves with the G curves of 
the system was used as a criterjion for the solution of the stability. 
equations, (7-1) and (7-2). 

Although this approach is' convenient in the case of a single 
nonlinearity, it becomes very /cumbersome, if not impossible to use, 
when more than one nonlinearity exists, except in a few special 
cases [8]. For example, witli two coupled LST axes, the stability 
equations in the continuous 1 case may be of the form. 


1 + G ] (jw)N(A) + G 2 (jw)N^(A) = 0 


Clearly, a graphical solution for Eq. (7-3) is impractical , and 


(7-3) 
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the situation very quickly deteriorates if more axes are added, or 
if the amplitudes of oscillations in the different axes are not the 
same. 

In view of this situation, it appears worthwhile to consider an 
alternate method for solving the stability equations of the type of 
(7-1) through (7-3). In this chapter the results obtained by solving 
Eqs. (7-1) and (7-2) for the single-axis LST system by means of an 
iterative numerical method are presented. This method has been 
successful and provides numerical solutions to Eqs. (7-1) and (7-2), 
these solutions being identical to those obtained by the earlier 
graphical method. The method is promising and can be directly extended 
to the more complicated cases, as in Eq. (7-3). 


7-2 The Numerical Method 

Consider the set of two nonlinear equations 
f-| (*i > *2) = ® 

f 2 (x r x 2 )=0 (7-4) 

or in more compact notation 

F(x) = 0 (7-5) 

where 


■ V 



' X 1 ' 



x = 


. f 2 



. X 2 . 


An algorithm for an iterative solution of Eqs. (7-4) or (7-5) is 

x k+1 = x k - CF x (x k )]"V(x k ) 


( 7 - 7 ) 
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k 

where is the value of x at the kth step, and £ x is the Jacobian 


r 3f 1 


ax^ 

ax 2 

af 2 

af 2 

ax 1 

ax^ 


The algorithm in Eq. (7-7) is the well-known Newton's method for 
the multivariable case, and can be shown to locally converge to the 
solution x* of Eq. (7-5) if 

i) £(x*) = 0 

and 

ii) E x (x*) exists. 

Recently, a numerical method has been proposed by Brown [9, 10] 
which is similar but computationally more efficient than Newton's 
method and still posesses the same convergence properties. This 
method is now used to solve Eqs. (7-1) and (7-2). 

7.3 The Continuous Case 

The stability equation in the continuous case is as in Eq. (7-1) , 

1 + N(A)G(ju>) = 0 (7-1) 

where 

oj is the frequency in rads/ sec 

A Is the amplitude of the input sinusoid to the nonlinearity 
N Is the describing function 

and G is the system transfer function seen by the nonlinearity. 
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Since N and G are both complex quantities, they may be written as 
N(A) = N r (A) + jNj(A) 

G(jo)) = G R (ju>) + jGj(ju) (7-9) 

where N R , Nj, G R , and Gj are all real quantities.. 

Substituting Eq. (7-9) into Eq. (7-1) yields 

1 + (N r + jNj)(G r + - jGj) = 0 (7-10) 

Collecting the real and imaginary terms in Eq. (7-10) and equating 
them to zero gives 

1 '+ N R (A)G R (joj) - N I .(A)G I (jw) = 0 

N r (A)G i (» + Nj (A)G R ( jo)) = 0 (7-11) 

Using the notation of Eq. (7-5), Eq. (7-11) becomes 



f 1 = T + N r G r - NjGj 
f 2 = N R G I + N I G R 


X| = 0) 

x 2 = A (7-13) 

Thus, Eq. (7-11) represents a system of two simultaneous nonlinear 
equations in two unknowns. This system is solved by the proposed numerical 
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method with the following system parameters 


t gfo = 0.1 ft-ib ; 


H = 600 ft-lb-sec 



2.1 ft-lb-sec 2 I 

i 

5758.35 

1371.02 

216 ft-lb/rad/sec 
8700 ft- lb/ rad 
10 5 ft-lb-sec 2 


\ 





With y = 1.38 x 10 7 , it is known from the graphical results [2] 
that two solutions to Eq. (7-11) exist. Figures 7-1 and 7-2 show the 
numerical iterations for several initial solutions. The two solutions 
are 


1 . w = 4.27 rad/sec 

A = 5.57 x io ' 6 rad 

2. a) = 1.748 rad/sec 

A = 4.45 x IQ' 7 rad 


Although solution 2 is an unstable equilibrium and solution 1 is 
a stable equilibrium, the numerical method does not differentiate between 
than. It converges on the solution in whose domain of attraction the 
inital solution is chosen. Figure 7-3 shows a graphical interpretation 
of the domains of attraction of the two solutions. If the initial 
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4 

0.4259855835D 01 

0.5464563603D-05 

5 

0 .42719361 1 ID 01 

0 . 557 1641 673D- 05 

6 

0.4272145453D 01 

0.557346 7 159D- 05 

7 

0 .42721 45393D 01 

0.5573465998D-05 

ND= 1 FREQUENCY® 4.27215D OORAD/SEC AMPLITUDE® 5.57347P-06 

-L 

— ^ " 


—Figure 7-1 

. Numerical iterations in the continuous case; 


Y = 1.38 x TO 7 . 
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LST SYSTEM-1 

NUMERICAL SOLUTION 

OF 1+N(A)GCW)=0. 

GAMMA* 1.38000D 07 



ITERATION 

FREQUENCY 


AMPLITUDE 

0 

0 .500 00 0000 OD 

01 

0.5000000000D-05 

1 

0.433 0324 789D 

01 

0 .571607 0923D- 05 

c. 

0.4272669199D 

01 

0. 5575178 057D-05 

3 

0.4272145539D 

01 

0.5573465283D-05 

4 

0.42721 45393D 

01 

0.55734660 OOD- 05 

NB= 1 FREQUENCY= 4.27215D 

OORAD/SEC AMPLITUDE* 5.57347D-06 

LST SYSTEM- 

NUMERICAL SOLUTION 

OF 1+N(R>G(W)=0. 

GAMMA* 1.38000D 07 



ITERATION 

FREQUENCY 


AMPLITUDE 

0 

0.5000000000D 

01 

0.1 OOOOOOOOOD-05 

1 

0. 398628 1834D 

01 

0.2781931 045D-05 

2 

0.408101 1589D 

01 

0.4234094904D-05 

3 

0. 42423334 OOD 

01 

0.5306434465D-05 

4 

0 .4270849057D 

01 

0.5562073786D-05 

5 

0.4272143796D 

01 

0.5573455629D-05 

6 

0 . 4272 1 45393D 

01 

0.5573466007D-05 

ND= 1 FREQUENCY* 4.272151 

i OORAD/SEC AMPLITUDE* 5.57347D-06 

LST SYSTEM- 

NUMERICAL SOLUTION 

OF 1+N<A>G(«)=0. 

GAMMA* 1 .38 GOOD 07 



ITERATION 

FREQUENCY 


AMPLITUDE 

0 

0.3000000000D 

01 

0.100000000 OD- 05 

1 

0.41 0676 1584D 

01 

0 .3045 086482D- 05 

2 

0.411 0062736D 

01 

0.4414687420D-05 

3 

0.4249871 8 06D 

01 

0 . 5373028 097D-05 

4 

0.4271417631D 

01 

0 . 55670843 1 5D- 05 

5 

0.4272145079D 

01 

0.5573465364D-05 

6 

0.4272145393D 

01 

0.5573465999D-05 

ND= 1 FREQUENCY* 4.272151 

l GORAD/SEC AMPLITUDE* 5.57347D-06 


LST SYSTEM- 

NUMERICAL SOLUTION OF 

1 ♦N < A ) G < W > * 0 . 

GAMMA* 1.380 OOD 07 



ITERATION 

FREQUENCY 


AMPLITUDE 

0 

0..2000000000D 

01 

0.5000000000D-06 

1 

0 . 1633524697D 

01 

0.28 1 297 0363D- 06 

c 

0.1 71293 0384D 

01 

0.4060843547D-06 

o 

0.1 745466 022D 

01 

0.4418355628D-06 

4 

0.1 748 094452D 

01 

0.4451 0988 05D- 06 

5 

0.1 748 1 04365D 

01 

0.44511 824 16D- 06 

ND= 1 FREQUENCY* 1.74810D 

OORAD 

/SEC AMPLITUDE* 4.451 18D- 07 


Figure 7-2. Numerical iterations in the continuous case; 
= 1.38 x 10 7 . 


Y 





Figure 7-3. Domains of attraction of the two solutions with 
Y = 1.38 x 10 7 . 
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conditions are chosen outside of the domain of attraction, the method 
does not converge. 

With changes in y, the solution points shift, and Fig. 7-4 shows 

Q O 

the stable solutions with y = 0.69 * 10 and 1.38 x 10 respectively. 

7.4 The Discrete Case 

In this case the input to the nonlinearity is assumed to be 
e(t) = Acos(a)t + <p) (7-14) 

and the stability equation is 

1 + N(A, <p, n)G(T, n, z) = 0 (7-15) 

where 

T is the sampling period of the system 

n is the order of oscillation, i.e., the period of oscillation is Tn. 

A is the amplitude of the sinusoidal input to the nonlinearity 

4> is the phase of this input relation to the sampling process 

z = exp(j2ir/nT) 

6 is the z-domain transfer function seen by the nonlinearity 

N is the discrete describing function of the nonlinearity i 

To maintain consistency, Eq. (7-15) can have two and only two 
variables; thus, n and 4 > are assumed to be fixed parameters, and A and 
T are the two variables. With each value of n and <)>, a solution of 
Eq. (7-15) is desired. 

As in the continuous case, define 

G(T) = G r (T) + jGj(T) 

N(A) = N r (A) + j N j (A) 


(7-16) 
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LST SYSTEM-NUMERICAL SOLUTION OF 1+N< A>G<W>=0 . 
GAMMA= 6.90000D 07 


ITERATION 

FREQUENCY 


AMPLITUDE 

0 

0.7000000000D 

01 

0.1 00 0000 000D- 04 

1 

0.5323793 06 OD 

01 

0.1 064495497D-04 

•-> 

0.5203573337D 

01 

0.81 13476180D-05 

3 

0.522 09 17925D 

01 

0.8556741622D-05 

4 

0 . 522 1 777385D 

01 

0.85S0953482D-05 

5 

0.5221778880D 

01 

0.858 0997 11 3D- 05 


ND= 1 FREQUENCY* 5.22178D OORAD'SEC AMPLITUDE* 8.581 00D-06 


LST SYSTEM-NUMERICAL SOLUTION OF 1«-N<A)G< U> = 0 . 
GAMMA* 6.90000D 07 


ITERATION FREQUENCY 

0 0.50000 0000 OD 01 

1 0.520480881 OD 01 

2 0 . 522 127061 9D 01 

3 0 .5221 778574D 01 

4 0.5221778879D 01 

ND= 1 


AMPL I TUDE 
0 . 1 000000000D-04 
0.8217060863D-05 
0 .8586884 03 ID- 05 
0.8580987575D-05 
0.8580997080D-05 

8.581 00D- 06 


FREQUENCY* 5.22178D OORAD/SEC AMPLITUDE* 


LST SYSTEM-NUMERICAL SOLUTION OF 1+N<:A>6<W) = 0 . 
GAMMA* 1 .380 GOD 08 


ITERATION 

FREQUENCY 


AMPLITUDE 

0 

0.7000000000B 

01 

0 . 1 OOOOOOOOOD-04 

1 

0. 5463694 124D 

01 

0 . 1 074650663D-04 

o 

c 

0.5412668486D 

01 

0.3894201 4 14D- 05 

3 

0. 54 1931 659 OD 

01 

0 .9155958622D- 05 

4 

0.541 9498774D 

01 

0 . 9163497 173D- 05 

KT 

0 .5419498719D 

01 

0.9163495937D-05 


HD* 1 FREQUENCY* 5. 41950D OORAB'SEC AMPLITUDE* 9.16350D-06 


LST SYSTEM-NUMERICAL SOLUTION OF 1 ♦N<A>G<U>=0 . 
GAMMA* 1.38000D 03 


ITERATION FREQUENCY AMPLITUDE 


0 

0.5 OOOOOOOOOD 

01 

0 . 1 OOOOOOOOODr 04 

1 

0.5396836346D 

01 

0. 899265382 OD- 05 

2 

0.541 944 1243D 

01 

0.91615031 19D- 05 


0 .5419498763D 

01 

0.9163497435D-05 

4 

0 .54 1 949371 9D 

01 

0 . 9 1 63495936D- 05 

FREQUENCY* 5.41950D 

0 ORAD 

■SEC AMPL I TUDE 


Figure 7-4. Numerical iterations in the continuous case; 
Y = 0.69 x 10 8 and 1.38 x io 8 . 
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Substituting Eq. (7-16) into Eq. (7-15) and separating the real and 
imaginary parts yields 

f l^ x l’ x 2^ 0 

F(x) = = (7-17) 

^2^ x i * x 2^ 0 

. 

where 

f l = 1 + N R G r - 

f 2 = N r Gj + NjG r 

X 1 = T 

x 2 « A (7-13) 

Equation (7-17) now represents two equations in two unknowns and 
can be solved by the proposed numerical method. With the same system 
parameters as in the continuous case, and y - 1.38 x 10 7 , a solution 
for n = 20, <J> = 0 is obtained. The initial solution in this case is 
obtained from a knowledge of the continuous system solution. This 
solution, n =20 and 4> = 0, is now used as an initial guess to determine 
the solution for n = 18, <j> = 0, which is then used as an initial solution 
to obtain the solution for n = 16, 4> = 0, and so on. Figures 7-5, 7-6 
and 7-7 show the iterations for n = 20 through n =4. , Due to the different 
characteristics of the odd-n and even-n solutions, the decrement or 
increment of n is made in multiples of 2. This way, the solution for 
all n can be obtained. 

Once a solution is available for a particular n at <j> = 0, the 
<t>-space can be spanned by slowly varying f and using the previous 
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LST SYSTEM-NUMERICAL SOLUTION OF DISCRETE DESCRIBING FUNCTION 
GAMMA® 1.38000D 07 N= 20 PHI= 0.0 

ITERATION SAMPLING PERIOD AMPLITUDE 

0 0.6000000000D-01 0 . 1 0000000000-04 

1 0.6635022584D-01 0 .71 68351 284D-03 

2 0 . 6684849597D- 0 1 0.739159141 ID-03 

3 0.6684184809D-01 0 .7402030884D-05 

4 0.6684185386D-01 0 .7402032739D-05 

ND® 1 SAMPLING PERIOD* 6 .68419D-02SEC AMPLITUDE* 7.40203D-06 


LST SYSTEM-NUMERICAL SOLUTION OF DISCRETE DESCRIBING FUNCTION. 
GAMMA* 1.38000D 07 N* 18 PHI® 0.0 

ITERATION SAMPLING PERIOD AMPLITUDE 

0 0.6684185385D-01 0 .7402032738D-05 

1 0 . 725898 1714D-01 0 .783 18 0082 0D- 05 

2 0 . 7332495463D- 0 1 0 .774 19651 86D- 05 

3 0.7333309288D-01 0 .774 0553676D-05 

4 0.7333308321D-01 0 .7740556323D-05 

ND* I SAMPL ING PER I OD* 7.3333 ID- 02SEC AMPL I TUDE* 7 . 74 056D- 06 


LST SYSTEM-NUMERICAL SOLUTION OF DISCRETE DESCRIBING FUNCTION 
GAMMA* 1.38000D 07 N= 16 PHI* 0.0 

ITERATION SAMPLING PERIDD AMPLITUDE 

0 0. 73333 08323D- 01 0 .774 03563 18D- 03 

1 0.8035141927D-01 0 .8308808703D-05 

2 0.813681 0323D- 01 0 .8182471 046D-05 

3 0.81331 85695D- 0 1 0 .8180174756D-05 

4 0 . 8 1 38 1 84229D- 0 1 0 .8180178901D-05 

ND* 1 SAMPLING PERIOD* 8. 1381 8D-02SEC AMPLITUDE* 8. 1801 8D- 06 


LST SYSTEM-NUMERICAL SOLUTION OF DISCRETE DESCRIBING FUNCTION 
GAMMA* 1.38 00 0D 07 N= 14 PHI* 0.0 

ITERATION SAMPLING PERIOD AMPLITUDE 

0 0 .8 1381 8423 ID- 01 0 .8180178895D-05 

1 0.9022050932D-01 0 .8961 126589D-05 

2 0 . 9 175 1 28977D- 01 0 .8758961 044D- 05 

3 0.91 7788 0263D- 0 1 0. 875466 0548D- 03 

4 0. 91 7787785 0D- 01 0 .875466 7872D- 05 • 

ND* 1 SAMPLING PERIOD* 9 . 17788D-02SEC AMPLITUDE* 8.75467D-06 


Figure 7-5. Numerical iterations in the discrete case; 
Y = 1.38 x 10 7 , $ =0, variable n. 
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LST SYSTEM-NUMERICAL SOLUTION OF DISCRETE DESCRIBING FUNCTION 
GAMMA* 1 .38000D 07 N= 12 PHI= 0.0 

ITERATION SAMPLING PERIOD AMPLITUDE 

0 0.9177877353D-01 0.87546678610-05 

1 0 . 1 033226304D 00 0.88918766750-05 

£ 0 . 1 058887009D 00 0 .3522086634D-05 

3 0. 1 05956962 0D 00 0 .9512652996D-05 

4 0 . 1 059569299D 00 0 .9512665872D-05 

ND- 1 SAMPLING PERIOD- 1 . 05957D-01SEC AMPL I TUDE= 9.51267D-06 


LST SYSTEM-NUMERICAL SOLUTION OF DISCRETE DESCRIBING FUNCTION 
GAMMA* 1 .38000D 07 N= 10 PHI= 0.0 

ITERATION SAMPLING PERIOD AMPLITUDE 

0 0.1 0595693 00D 00 0 .9512665855D-05 

1 0 . 121601231 ID 00 0 . 1 1 31 7 08575D-04 

2 0. 126490871 ID 00 0 . 1 056350696D- 04 

3 0. 1267133857D 00 0 . 1 053796573D-04 

4 0 . 1267135528D 00 0. 1 053796318D-04 

ND= 1 SAMPLING PERIOD- 1 .26714D-01SEC AMPLITUDE- 1.05380D-05 


LST SYSTEM-NUMERICAL SOLUTION OF DISCRETE DESCRIBING FUNCTION 
GAMMA- 1.38000D 07 N= 3 PHI- 0.0 

ITERATION SAMPLING PERIOD AMPLITUDE 

0 0 . 1 2671 35526D 00 0 . 1 05379632 0D- 04 

1 0 . 1484753967D 00 0 . 1378146933D-04 

2 0 . 1 59339 1474D 00 0 . 121 0637147D-04 

3 0 . 1604043798D 00 0 . 1200135682D-04 

4 0.16041 09954D 00 0 . 1200077001D-04 

ND- 1 SAMPLING PERIOD* 1 .6041 1D-01SEC AMPLITUDE- 1 .20008D-05 


LST SYSTEM-NUMERICAL SOLUTION OF DISCRETE DESCRIBING FUNCTION 
GAMMA* 1.38000D 07 N= 6 PHI- 0.0 

ITERATION SAMPLING PERIOD AMPLITUDE 

0 0 . 16041 09885D 00 0 . 12000771 06D-04 

1 0.1897914063D 00 O. 1878650442D-04 

2 0.2182241830D 00 0 . 1 497205004D-04 

3 0 .2272750301 D 00 0 . 141 0272682D-04 

4 0 . 227737 1225D 00 0 . 1408970958D-04 

5 0.2277376514D 00 0 . 1408969070D-04 

ND- 1 SAMPLING PERIOD- 2 .27738D-01SEC AMPLITUDE- 1 .40897P-05 


Figure 7-6. Numerical iterations in the discrete case; 
Y = 1.38 x 10 7 , <p = 0, variable n. 
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LST SYSTEM-NUMERICAL SDLUTION OF DISCRETE DESCRIBING FUNCTION 
GAMMA* 1.38 00 OD 07 N* 4 PHI= 0.0 

ITERATION SAMPLING PERIOD flNPLITUDE 


0 

0.227737 6508D 

00 

0 ; 1 4 08969 073D- 04 

1 

0.2543144363D 

00 

0 . 2792562588D- 04 

2 

0 .31 19051493D 

00 

0.258422275 ID- 04 

3 

0.41 8637 0592D 

00 

0 . 1482765201D-04 

4 

0.5071 624 02 OD 

00 

0.15571921 98D- 04 

5 

0 . 5354877655D 

00 

0.1 58493 1294D- 04 

6 

0. 537 14691 61 D 

00 

0 . 15870703930-04 

7 

0.53715 03 033D 

.00 

0.1 587 076 01 3D- 04 


ND= 1 SAMPLING PERIOD- 5. 371 SOD- 01 SEC AMPLITUDE- 1 .587 08D- 03 


Figure 7-7. Numerical iterations in the discrete case; 
Y = 1-38 x TO 7 , 4> = 0, variable n. 
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solution as the initial guess for the next solution. The range of <f> 
is 2ir/n for n even and ir/n for n odd. Figures 7-8 through 7-11 show 
the iterations for n = 3 and variable <j>. A plot of all solutions of 
T and A for n = 3 and various <f>'s is shown in Fig. 7-12. Here the 
plots for n = 4, 6, 8 and 10 are also shown. As n increases the range 
of <J> decreases and the plots of T and A with $ as a parameter shrink 
to single points. Figure 7-13 shows the plot of T and A with <t> = 0 
for n = 10 through n = 100. 

The plots of Figs. 7-12 and 7-13 can be used to determine the 
frequency (2ir/nT) and amplitude (A) of self-sustained oscillations when 
the sampling period T is given. They provide very concise information 
on the conditions under which self -sustained oscillations can occur in 
the discrete-data system. The corresponding solutions by the graphical 
method require overlapping the -1/N and G curves for each n and checking 
for intersection or containment. 

As n increases beyond 100, the discrete solution asymptotes 
towards the continuous solution. 

7.5 Conclusions and Extensions 

The results of the previous sections have demonstrated the ease 
and effectiveness of the proposed numerical method. It provides the 
exact frequency and amplitude of oscillation in comparison to the 
approximate values obtained by the graphical technique. In the discrete 
case particularly, the numerical approach provides the information in 
a more convenient form. The frequency and amplitude combinations which 
can exist for each value of sampling period are available from a single 


curve. 




The major limitation of the numerical method is the convergence 
characteristics of the solutions. Unless an adequate initial solution 
is available, no useful information can be obtained. The initial 
solution must be near to or within the domain of attraction of the 
exact solution. Once a solution is obtained, it is simple to slowly 
vary the parameters and obtain all the desired solutions. 

The method is unrestricted to the form of the stability equation 
or the number of nonlinearities present. If adequate initial solutions 
can be selected, this approach can yield useful results with. the more 
complicated multi pie -axis models of the LST. 
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LST SYSTEM-NUMERICAL SOLUTION OF DISCRETE DESCRIBING FUNCTION 
GAMMA* 1.38000D 07 N= 3 PHI= 6.00000D 01 
ITERATION SAMPLING PERIOD AMPLITUDE’ 

0 0.1000000000D 01 0.2000000000D-04 

1 0.1206736268D 01 0 .2180807968D-04 

2 0.1246507649D 01 0 .22 0096637 0D- 04 

3 0 . 1247557574D 01 0 .2201 127396D-04 

4 . 0 . 1247557444D 01 0 .2201 127250D-04 

ND= 1 SAMPLING PERIOD* 1 .24756D OOSEC AMPLITUDE* 2.20U3D-03 


LST SYSTEM-NUMERICAL SOLUTION OF DISCRETE DESCRIBING FUNCTION 
GAMMA* 1.38000D 07 N= 3 PHI* 5.00000D 01 : 

ITERATION SAMPLING PERIOD AMPLITUDE 

0 0 . 1247557444D 01 0 .2201 127230D-04 

1 0.7576140716D 00 0 .221 0012199D-04 

2 0 . 8593652656D 00 0 .221201 51 79D-04 

3 0.871 0667 172D 00 0 .2219954845D-04 

4 0.871 1804237D 00 0 .2220037576D-04 

5 0.871 1803498D 00 0 .2220037467D-04 

ND= 1 SAMPLING PERIOD* 8. 71 180D-01SEC AMPLITUDE* 2.220 04D- 05 


LST SYSTEM-NUMERICAL SOLUTION OF DISCRETE DESCRIBING FUNCTION 
GAMMA* 1.38000D 07 N= 3 PHI* 4.00000D 01 
ITERATION SAMPLING PERIOD AMPLITUDE 

0 0.871 1803498D 00 0 .2220037468D-04 

1 0.6048096982D 00 0 .22991 03460D-04 

£ 0 .650 0626463D 00 0 .23343301 00D- 04 

3 0. 652778725 OD 00 0 .234266340 ID- 04 

4 0 . 65278444 1 6D 00 0 .234271 0075D-04 

ND* 1 SAMPLING PERIOD* 6. 52784D-01SEC AMPLITUDE* 2. 3427 ID- 05 


Figure 7-8. Numerical iterations in the discrete case; 
y = 1.38 x io 7 , n = 3, variable 4>. 
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LST SYSTEM-NUMERICAL SOLUTION OF DISCRETE DESCRIBING FUNCTION 
GAMMA* 1.38000D 0? N= 3 PHI= 3.00000D 01 
ITERATION SAMPLING PERIOD AMPLITUDE 

0 0.6527844384D 00 0 .£34271 0022D-04 

1 0.4972491700D 00 0 .2522386 042D- 04 

£ 0.51 85445957D 00 0 .2594052 161D-04 

3 0 .5 191 962437D 00 0 .2604023247D-04 

4 0.5191954834D 00 0 .2604065065D-04 

ND= 1 SAMPLING PERIOD* 5 . 19195Dr01SEC AMPLITUDE* 2.60407D-05 


LST SYSTEM-NUMERICAL SOLUTION OF DISCRETE DESCRIBING FUNCTION 
GAMMA* 1.38000D 07 N* 3 PHI* 2.00000D 01 
ITERATION SAMPLING PERIOD AMPLITUDE 

0 0.5191954849D 00 0 .2604065029D-04 

1 0.4479887233D 00 0 .27761 94476D- 04 

2 0.4542599663D 00 . 0 .2331780647D-04 

3 0.454254854 ID 00 0 .2835057365D-04 

4 0.454254840 ID 00 0 .2835057996D-04 

ND* 1 SAMPLING PERIOD* 4 .54255D-01SEC AMPLITUDE* 2.83506D-05 


LST SYSTEM-NUMERICAL SOLUTION OF DISCRETE DESCRIBING FUNCTION 
GAMMA* 1.38000D 07 N* 3 PHI* 1.00000D 01 
ITERATION SAMPLING PERIOD AMPLITUDE 

0 0. 45425484 02D 00 0 .2835057996D-04 

1 0.6655 15 049 ID 00 0 . 1 157357946D-04 

2 0 . 73038895 06D 00 0 . 16422 1382 ID- 04 

3 0.6664793457D 00 0 . 1 959748765D-04 

4 0.6559819499D 00 0 .2045437953D-04 

5 0. 655378 0327D 00 0 .2049946525D-04 

6 0 . 65537775 02D 00 0 .2049953882D-04 

ND= 1 SAMPLING PERIOD* 6 .55378D-01SEC AMPLITUDE* 2. 04995D- 05 


Figure 7-9. Numerical iterations in the discrete case; 
Y = 1.38 x 10^, n = 3, variable (j>. 



2 


3 


4 


5 


Amplitude (*lCf^) 

Figure 7-12. Amplitude and frequency (2Tr/nT) of self sustained 
oscillations for various sampling periods in the 
discrete case; y = 1 . 38 * 10^ , n = 3 through n = 10. 
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L ST SYSTEM-NUMERICAL SOLUTION OF DISCRETE DESCRIBING FUNCTION 
GAMMA® 1.38000D 07 N= 3 PHI® 4.00000D 00 
ITERATION SAMPLING PERIOD . AMPLITUDE 

0 0. 944392 1729D 00 0 .200043 1869D- 04 

1 0 . 1 098499702D 01 0 .2069887736D-04 

2 0 . 1 120422983D 01 0 .2075227791D-04 

3 0. 1 120755444D 01 0 .2075252072D-04 

4 0. 1 120755263D 01 0.2075252039D-04 

ND= 1 SAMPLING PERIOD= 1.12076D OOSEC AMPLITUDE® 2.07525D-05 


LST SYSTEM-NUMERICAL SOLUTION OF DISCRETE DESCRIBING FUNCTION 
GAMMA® 1.38000D 07 N= 3 PHI= 2.00000D 00 
ITERATION SAMPLING PERIOD AMPLITUDE 

0 0 . 1 120755263D 01 0 .2075252 039D- 04 

1 0 . 1243920308D 01 0 .21595141 12D-04 

2 0 . 1255495749D 01 0 .2162959060D-04 

3 0. 1255574083D 01 0 .2162960044D-04 

ND® 1 SAMPLING PERIOD® 1 .25557D OOSEC AMPLITUDE® 2.16296D-05 


Figure 7-11. Numerical iterations in the discrete case; 
Y = 1.38 x 10^, n = 3, variable <J>. 
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LST SYSTEM-NUMERICAL SOLUTION 
GAMMA* 1.38000D 0? N= 3 
ITERATION SAMPLING PERIOD 

0 0. 600000000 OD' 00 

1 0.653011 1015D 00 

£ 0.6553859033D 00 

3 0. 655377732 OD 00 

ND= 1 SAMPLING PERIOD* 6.553 


OF DISCRETE DESCRIBING FUNCTION 
PHI* 1.00000D 01 
AMPLITUDE 
0.2000000000D-04 
0.2038544865D-04 
0 . £ 049873663D- 04 
0.2049953953D-04 

'8D-01SEC AMPLITUDE* 2.04995D-05 


LST SYSTEM-NUMERICAL SOLUTION 
GAMMA* 1.38 00 OD 07 N= 3 
ITERATION SAMPLING PERIOD 

0 0.6553777518D 00 

1 0.77 03562 04 ID 0 0 s 

£ 0.7841 108478D 00 

3 0.784207581 ID 00 

4 0.7842075975D 00 


OF DISCRETE DESCRIBING FUNCTION 
PHI* 8.00000D 00 
\ \AMPLITUDE 
0 .2049953873D-04 
^0.1 958935672D- 04 
0.1 98£ 12191 9D- 04 
0 . 1 982549843D- 04 
0.1 982549437D- 04 


ND= 1 


SAMPLING PERIOD* 7 .842 OSD- 01 SEC 


ArfpL\lTUDE= 




1 . 98255D- 05 


LST SYSTEM-NUMERICAL SOLUTION 
GAMMA* 1 .38 0.0 OD 07 N= 3 
ITERATION SAMPLING PERIOD 

0 0 . 7842075974D 00 

1 0. 923581 0128D 00 

£ 0.944071261 ID 00 

3 0. 9443923 074D 00 

4 0. 944392 1728D 00 


OF DISCRETE DESCRIBING FUNCTION 
PHI* 6.00000D 00 
AMPLITUDE 
0 . 1 982549438D-04 
0 . 1990482637D-04 
0.2000313242D-04 
0.2000432019D-04 
0.2000431 868D- 04 


ND= 1 SAMPLING PERIOD* 9 .44392D-01 SEC AMPLITUDE* 2.00043D-05 


Figure 7-10. . Numerical iterations in the discrete case; 
Y = 1.38 x io 7 , n = 3, variable <t>. 
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Figure 7-13. Amplitude and frequency (2ir/nT) of self sustained 
oscillations for various sampling periods in the 
discrete case; y = 1.38 x 10 7 , n = 10 through n = 100. 
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